SINGLE PARTICLE SLOW DYNAMICS OF CONFINED WATER 
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Molecular dynamics simulations of SPC/E water confined in a Silica pore are presented. The 
simulations have been performed at different hydration levels and temperatures to study the single- 
particle dynamics. Due to the confinement and to the presence of a hydrophilic surface, the dynamic 
behaviour of the liquid appears to be strongly dependent on the hydration level. On lowering temper- 
ature and/or hydration level the intermediate scattering function displays a double-step relaxation 
behaviour whose long time tail is strongly non-exponential. At higher hydrations two quite distinct 
subsets of water molecules are detectable. Those belonging to the first two layers close to the sub- 
strate suffer a severe slowing down already at ambient temperature. While the behaviour of the 
remaining ones is more resemblant to that of supercooled bulk SPC/E water. At lower hydrations 
and/or temperatures the onset of a slow dynamics due to the cage effect and a scenario typical of su- 
percooled liquids approaching the kinetic glass transition is observed. Moreover, for low hydrations 
and/or temperatures, the intermediate scattering function clearly displays an overshoot, which can 
be assigned to the so called "Boson Peak" . 



I. INTRODUCTION 

Water plays a most fundamental role on earth and its 
anomalous properties as a function of pressure and tem- 
perature are the subject of a longstanding scientific de- 
bate. Nevertheless there are still many crucial questions 
that remain to be answered. Some of the most impor- 
tant are described throughout all the present issue that 
is dedicated to metastable water Q|. In particular, the 
supercooled region has attracted a wide attention in the 
last few years. Thermodynamic as well as transport, dy- 
namic, properties are not completely understood in this 
region. The main reason is that experimentally water can 
be supercooled only down to Tjj ~ 236K. Then, because 
of homogeneous nucleation due to impurities, the liquid 
is driven toward the crystal phase ||. Therefore studies 
of the liquid both through classic computer Molecular 
Dynamics, MD, and through theoretical models offer us 
ways to penetrate deep inside that region in the attempt 
to release the information that is not directly accessible 
to experiments. 

In particular the normal limit of supercooling and the 
possibility of vitrification of water is a fundamental point 
that necessitates to be clarified. A recent experiment sus- 
tains the hypothesis that the amorphous phase can be 
connected to the normal liquid phase through a reversible 
thermodynamic path Q| . MD simulations of bulk super- 
cooled SPC/E jj] water showed a kinetic glass transition 
as predicted by mode coupling theory (MCT) || at a 
critical temperature Tq ~ Ts ||, where T5 is the singu- 
lar temperature of water Q], which is Ts = 228K or, for 
SPC/E, 49 degrees below the temperature of maximum 
density. 

Within this framework a comparison of the behavior 
of the bulk liquid with the same liquid in a confined 
environment is highly interesting since both the phase 



diagram and the dynamic behaviour of confined water 
could present an analogy with the same liquid in the 
bulk phase. The outcomings of this kind of researches 
can have therefore important experimental implication 
for the bulk. In particular experimentally forbidden re- 
gions of the phase diagram of bulk water could become 
accessible through the study of confined water. On the 
other hand a study of the modification of confined water 
with respect to the bulk is highly interesting for the de- 
velopment of both biological and industrial applications. 

It has been inferred from experiments, although not 
conclusively proved, that confining water could be equiv- 
alent to the supercooling of the bulk . In this paper 
I will inquire on how far this analogy can be pushed for 
SPC/E water confined in a silica pore and its bulk phase. 

Among the different systems studied experimentally 
water confined in porous Vycor glass is one of the most 
interesting with relevance to catalytic processes and en- 
zymatic activity. Vycor offers in fact to water a well 
characterized network of cilidrical pores, is composed of 
simple molecules SiC>2, has a quite well characterized 
structure with a sharp distribution of pore sizes with an 
average diameter of ~ 40±5A, i.e. the same order of mag- 
nitude of many biological confining environments, and a 
strongly hydrophilic inner surface. For these reasons sev- 
eral experiments on water-in- Vycor have been performed 
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I shall show in the following the results concerning the 
single particle dynamics obtained from a series of MD 
simulations of SPC/E water confined in a cylindrical sil- 
ica cavity The pore has been modeled to rep- 
resent the average properties of Vycor pores. The simu- 
lated system from one side represents therefore a rather 
general confining environment and from the other offers 
the possibility of a direct comparison with experiments. 
The results shown in the following are focused on the role 
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of hydration level on slow dynamics of confined water in 
comparison with the bulk properties upon supercooling. 

In the next section the details of the model of the pore 
are described together with the simulation details. In the 
third section the single particle dynamics of confined wa- 
ter is discussed in comparison with that of bulk SPC/E 
water. The last section concludes the paper with a dis- 
cussion on the results. 



II. MODEL AND SIMULATION DETAILS 



A single pore has been constructed for the simulations 
since in the range of timescales involved in this kind of 
dynamics water does not leak out of the pore. A cubic 
cell of vitreous S1O2 was built by melting a /3-cristobalite 
single crystal and then by quenching the system down to 
ambient temperature j[(|[l7|. The side of the cube is 
d ~ 70 A. Inside the cell a cylindrical cavity, the pore, 
with diameter of 40 A was created by removing all the 
atoms lying within a distance of R = 20 A from the z-axis 
passing through the center of the cube. The inner surface 
of the pore was then "corrugated" by removing from it 
all the Si atoms which were bonded to less than four O. 
At the surface oxygen atoms can be classified as bridging 
oxygens (bO) if they are bonded to at least two silicon 
atoms, non-bridging oxygens (nbO) otherwise. Then the 
dangling bonds of the nbO were saturated by attaching 
an hydrogen to each nbO. This procedure mimics the ex- 
perimental preparation of the sample of water-in- Vycor 
in which internal surfaces of dessicated Vycor are hy- 
drogenated before hydration. At the end of this process 
the cell contains 6400 silica atoms, 12500 bO and 230 
nbO. This last number yields a surface density of acidic 
hydrogen of 2.5 nm~ 2 in good agreement with the 2.3 
nrnT 2 obtained in the experiments |]l8| . The simulations 
have been performed in the NVE ensemble introducing 
water molecules into the pore. Water molecules inter- 
act among themselves via a SPC/E potential, and with 
Vycor atoms via an effective pair potential Jlq| . Both 
potentials have a Coulombic part plus a Lennard Jones 
term between oxygens. The parameters of the simulation 
are reported elsewhere [[l2]jlj| . The interactions cutoff is 
9 A and the shifted force method has been used. For 
the Coulomb interaction it was checked that with re- 
spect to the correct use of the Ewald summation this 
technique was not producing qualitatively different re- 
sults The timestep of the simulations was 2.5 fs. 
Equilibration, reached through a Berendsen thermostat 
pl| , was monitored via the time dependence of the po- 
tential energy and the total energy. Along the pore axis 
(z-direction) periodic boundary conditions are applied to 
the water molecules. Along the remaining directions (xy) 
water is confined. The glass is rigid. 
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FIG. 1. Snapshots from simulations of SPC/E water in- 
side a cylindrical cavity of diameter 40 A at ambient tem- 
perature. Only oxygens are shown, pw is the global water 
density in gj cm 3 (see Table |) . 



The molecular dynamics calculations have been per- 
formed for different numbers of water molecules, corre- 
sponding to different levels of hydration of the pore. The 
definition of hydration level of the pore deserves some 
consideration. Experimentally Vycor glass absorbs water 
up to 25% of its dry weight, this is defined as equilibrium 
or full hydration (hf ~ 0.25 g of water/ g of Vycor). The 
experimental density of water when confined in such a 
pore at full hydration is estimated to be 11% less than 
its bulk value at ambient conditions i.e. pw — 0.0297 
molecules/A 3 = 0.8877g/cm 3 @. A partially hydrated 
sample is then experimentally obtained by absorption of 
water in the vapor phase until the desired level of hy- 
dration is reached. For the designed pore of this sim- 
ulation the density of the full hydration is obtained for 
N\v = 2661, where N\y is the number of water molecules 
introduced in the pore. The five hydration levels investi- 
gated in the present work are reported in Table |. I will 
discuss here the results obtained for two temperatures: 
T = 298 K and T = 240 K. 

TABLE I. Hydration levels of the pore. Nw is the number 
of water molecules and pw the corresponding global density. 



(a) : based on estimated value for full hydration Nw = 2661 
molecules (see text). 
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500 


19% 


0.1687 


1000 


38% 


0.3373 


1500 


56% 


0.4971 


2000 


75% 


0.6658 


2600 


99% 


0.8788 
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FIG. 2. Radial density profiles normalized to the bulk 
for oxygen atoms at ambient temperature for the five hy- 
dration levels investigated: first curve on the top corre- 
sponds to pw = 0.8788(?/cm 3 , last on the bottom to 
pw = 0.1687gr/cm 3 . The vertical line schematically indicates 
the pore surface. Right inside the pore the density is different 
from zero because few molecules are trapped in the pockets 
of the corrugated pore surface (see text). 

The snapshots, see Fig.l, of the system at ambient 
temperature show that the simulated surface is strongly 
hydrophilic. At all the degree of hydration water ad- 
sorbs on the Vycor surface. A wide variety of water clus- 
ter is visible for the lowest global density investigated 
pw — 0.1687 g/cm 3 , which is lower than the estimated 
monolayer coverage (pw ~ 0.2219 g/cm 3 ). 

The density profiles of the oxygen atoms of confined 
water are shown in Fig. 2. It is observed already at lower 
hydrations the presence of a layer of water molecules wet- 
ting the substrate surface. At nearly full hydration two 
layers of water with higher than bulk density are evident. 
Few molecules are trapped inside small pockets close to 
the surface, which are a byproduct of the "sample prepa- 
ration process" . There is a strong tendency of water 
molecules close to the surface to form hydrogen bonds 
(HB) with the atoms of the substrate, in particular the 
hydrogens of water molecules with bO. As a consequence 
the HB network of water results to be strongly distorted 
close to the Vycor surface [ ^Ul^j . This is compatible 
with the findings of recent experiments on water-in- Vycor 
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III. SINGLE PARTICLE DYNAMICS OF 
CONFINED WATER 

Bulk SPC/E water, like many glass forming liquid, 
when supercooled, or when the pressure is increased, de- 
velops a diversification of relaxation times scales, one 



fast and one slow, well described by MCT ||. In the 
MCT H scenario a liquid approaches the glass transi- 
tion point with a dynamic behaviour mastered by the so 
called "cage effect" . The molecule or atom is trapped by 
the transient cage formed by its nearest neighbour, nn. 
After an initial ballistic regime it starts feeling the po- 
tential barrier of its nn. In this intermediate time region 
the particle is rattling in the cage formed by this poten- 
tial barrier. When the cage relaxes the particle is free to 
move and enters the normal diffusive regime. But as the 
liquid gets closer to Tc the cage relaxation time becomes 
longer and eventually approaches infinity in the idealized 
MCT. Below Tc cages are frozen and only hopping pro- 
cesses can restore ergodicity allowing the particle to move 
(extended MCT) |22]. It is important here to stress that 
the cage effect is usually due in liquids to the increase 
of density, i.e. of nn, and in this respect water plays a 
peculiar role since there is no substantial increase in the 
number of nn on supercooling. It is the increased stiff- 
ness of the HB network at low temperature that allows 
for a "cage" to be formed also in bulk supercooled water. 
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FIG. 3. MSD at ambient temperature of oxygens in the 
non-confined z-direction as a function of decreasing hydration 
from top to bottom. In the inset the diffusion coefficients 
are shown. They correspond to the various hydration lev- 
els (circles) and to bulk SPC/E water at ambient conditions 
(diamond). 

Signatures of a behaviour a la Mode Coupling can be 
found in all the correlators that have an overlap with the 
density. In the following I will analyze for confined water 
in particular the mean square displacement, MSD, and 
the single particle density-density correlation function in 
the Q,t space also called the intermediate scattering func- 
tion, ISF, or Fs(Q,t). The signatures of this diversifi- 
cation of relaxation times and therefore the shouldering 
of the relaxation laws are more evident for the Q values 
close to the peak of the structure factor that for oxygen at 
ambient temperature is Q = 2.25 A -1 . In particular the 
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long time tail of the ISF is predicted to have a stretched 
exponential behaviour when the system approaches the 
glass transition. 




t (PS) 

FIG. 4. ISF of oxygens at the peak of the structure fac- 
tor in the confined xy direction. The highly non-exponential 
long time tail could be fitted (dashed line) to eqjl| only for 
the lower hydration. 

In Fig. 3 the mean square displacement of oxygens is 
displayed in the non-confined z-direction at ambient tem- 
perature as a function of hydration level. After the initial 
ballistic diffusion, before entering the diffusive regime, a 
flattening of the curve at intermediate times is observed 
as the hydration level is decreased. In the inset the bulk 
diffusion coefficient and the confined ones are displayed. 
There is a substantial decrease of average mobility as the 
hydration level in the pore is lowered. Correspondingly 
the oxygens ISF at the oxygen-oxygen peak of the struc- 
ture factor, Fig. 4, displays a shouldering of the relaxation 
laws upon decreasing hydration level. All the tails of the 
correlators of Fig. 4 are highly non-exponential. None- 
the-less these ISF could not be fitted to the same formula 
used for bulk supercooled water 

F S (Q, t) = [1 - A(Q)] e~^ 2 + A(Q)e-^ (1) 

except for the lowest hydration. In eqjl] A(Q) — e~ a ® / 3 
is the Lamb-Mossbauer factor (the analogous of the De- 
bye Waller Factor for the single particle) arising from 
the cage effect, and r s and t; are, respectively, the short 
and the long relaxation times. For the fit of the ISF 
of pw = 0.1687 g/cm 3 shown in Fig. 4 a cage radius 
a ~ 0.44 A is obtained, which is similar to the radius 
obtained for bulk supercooled water where a ~ 0.5 A. 
The lower radius obtained for confined water may be due 
to the slightly higher density of water close to the sur- 
face, see Fig. 2. The short relaxation time t s ~ 0.14 ps is 
again comparable to the bulk value which is r s ~ 0.2 ps. 



ti ~ 356 ps and (3 = 0.35. The (3 associated with t\ 
is very low compared to the typical values obtained for 
the supercooled bulk. Moreover the lowest hydration 
correspond to less than the monolayer coverage so that 
the dynamics of this system is that of clusters of water 
molecules attached to the pore surface. For the lowest 
hydration it is also observed a bump around Q.lps that 
could be possibly related to the existence of the Boson 
Peak feature in the S(q,(j) This point will be dis- 

cussed at the end of the paragraph. For the remaining 
correlators of Fig. 4 no analytic function was found that 
could fit the strongly non-exponential tails. Due to the 
strong hydrophilicity of the pore a diversification of dy- 
namic behaviour is to be expected as we proceed from 
the pore surface to the center of the pore. In Fig. 5 the 
Fs(Q,t) for the highest global density at ambient tem- 
perature is split into the contribution coming from the 
two layers of water molecules closer to the pore surface, 
outer shells, and into the contribution coming from all 
the remaining ones, inner shells. The first two layers are 
defined according to the density profile, Fig. 2. The inner 
shell contribution could be perfectly fit to eq.[j] as shown 
in the figure, while the outer shells one decays to zero 
over a much longer timescale so that water molecules 
there behave already as a glass. From the fit it is ex- 
tracted (3 = 0.74, n = 0.75 ps t s = 0.21 ps. The fit 
shows a remarkable agreement and both the r values are 
similar to those of bulk water, while for bulk water at 
ambient temperature (3 = 1. This kind of analysis has 
been extended for a different hydration and a different 
temperature. 
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FIG. 5. ISF for oxygens at almost full hydration in the 
confined xy direction. Shell Analysis. The top curve is the 
contribution to the total ISF coming from the molecules in 
the two shells closest to the substrate, see Fig. 2. The bottom 
curve is the contribution from the molecules in the remaining 
shells. The central one is the total ISF. The dashed line is 
the fit to eqQ 
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t (ps) FIG. 7. t i (main picture) and /3 (inset) values extracted 

FIG. 6. Shell analysis for the supercooled regime. The from the fit of the lon S time tai j of the ISF to e lfl The con " 

dashed line is again the fit to eq. @. Note that the bump tmuous lme has slo P e 2 - The Q behaviour has been found m 

related to the BP is now evident only for the inner shell con- S lass formers close to the S lass transition, 
tribution. 



In Fig. 6 the ISF shell analysis is shown at the peak 
of the structure factor for T=240 K and p w = 0.4971 
g/cm 3 . Also in this case the fit to the stretched expo- 
nential of the tail of the inner shells contribution is re- 
markable. From the fit (3 = 0.62, 77 = 11 ps t s = 0.16 
ps are extracted. Again the 77 is similar to that of bulk 
water at the same temperature while the (3 here is much 
lower. In Fig. 7 the 77 and the (3 values extracted from 
the fits to eq.[l] done as a function of Q for T=240 K and 
pw — 0.4971 g/cm 3 are shown. The (3 value reaches a 
plateau value and the 77 values show a Q 2 dependence. 
Both these behaviours are typical of a glass former un- 
dergoing a kinetic glass transition and in particular the 
Q 2 behaviour has been observed for example in glycerol 
close to the glass transition |^3j . 

Let us now come back to Fig. 4 where for the lowest hy- 
dration a bump around 0.7 ps is clearly visible. Deeply 
supercooled bulk water displayed the same peak around 
0.35 ps HH^]. This overshoot has been observed in simu- 
lations of supercooled strong glass formers like Si02 [p3| . 
It has been attributed in literature to the existence of a 
Boson Peak (BP) in the frequency domain. The BP is 
an excess of vibrational modes present in many glasses. 
When this glassy anomaly appears in the liquid state it 
is considered as precursor of the glass transition. This 
peak in the ISF has been also alternatively related to a 
disturbance propagating in a finite box. If there are pe- 
riodic boundary conditions imposed (z direction for the 
present system) the disturbance would reenter the box 
after t = L/v s where v s is the sound velocity pq |. In 
spite of the possible existence of this finite size effect it 
has been recently proven that the BP can be detected by 
MD ||. 



Note that in this system the same bump is present also 
in the supercooled confined simulation, see Fig. 6, and it 
is more evident for the inner shells. So that the presence 
of a disordered surface from one side seems to enhance 
the temperature at which the BP shows up and also its 
intensity with respect to the bulk. From the other side 
since the silica glass is rigid in our simulation the water 
molecules in contact with the glass cannot vibrate much 
so that practically only the inner shells can sustain vi- 
brations. In Fig. 8 same correlators of Fig. 6 in the xy and 
z direction are shown. It is important to stress here that 
for this system the finite size effect in the xy directions 
is a feature also of the real pore. 




t(ps) 



5 



FIG. 8. Shell analysis for the supercooled regime in the 
confined (xy), non confined (z) and in an average random di- 
rection (central curve in each group of three). In particular: 
ISF from the two shells closest to the substrate (top curves), 
inner shells (bottom curves) and total (central curves). 

This bump is therefore likely to be observed for water- 
in- Vycor. Both for the t corresponding to the minima 
and the maxima of the bump appearing in the lowest 
curves of Fig. 8 the ratio of the t z /t xy is 1.44 against 
the 1.75 of the ratio between the lengths of the pore 
L z = 70A and L xy — AOA. Moreover, although less ev- 
ident, a BP can be observed also for the outer shells in 
Fig. 8 and the location is completely different with respect 
to that of the inner shells. So that in our ISFs there is 
no direct quantitative connection with a finite size effect 
that can be detected from the shape and the location of 
these bumps. It can be inferred from these data that 
both the size and the geometry and the presence of sur- 
face strongly influence the nature of this peak. Size and 
geometry signatures can be seen in the change of location 
of the peak in the ISF in xy and z directions, see Fig. 8. 
Signatures of a presence of a hydrophilic rigid surface can 
be related to the fact that the bump is more evident for 
the inner shells and is present already at room temper- 
ature. The changes of location and shape with respect 
to the bulk are therefore also probably influenced by all 
these causes. None-the-less, no quantitative statements 
can be made at this stage. Therefore both a frequency 
domain analysis and a instantaneous normal mode anal- 
ysis can help to shed light on this issue and they are 
currently in progress on this system. 



rable at the same temperatures but the /3 values of the 
confined water are much lower than those of the bulk 
always at the same temperature. About the so called 
Kolhraush exponent (3 can be said that in general it is 
different for different systems but apparently no special 
significance can be attributed to its numerical value [0. 
It is also important to mention here that SPC/E po- 
tential is considered one of the best potential for water 
upon supercooling, but these features are not produced 
by all potentials. For example the ST2 shows a "jump 
diffusion" behaviour upon supercooling p3 |. These two 
potential in particular are known to sandwich the be- 
haviour of experimental water so that ST2 water is less 
and SPC/E water is more structured than real water. 
On the other hand experimental signatures of MCT be- 
haviour have been found both in the bulk [p9| and in 
water-in- Vycor |^^] and also experimental signatures of 
a possible existence of a BP |}0| in water-in- Vycor have 
been detected. These results are encouraging and more 
simulations on this pore as a function of temperature 
and hydration are in progress for a full MCT test and a 
complete comparison with the bulk. 



V. ACKNOWLEDGEMENTS 

I wish to thank Mauro Rovere, Maria Antonietta Ricci 
and Eckhard Spohr for their contributions to this work, 
C. Austen Angell, Burkhard Ceil, Alberto Robledo and 
Francesco Sciortino for stimulating discussions on topics 
related to this paper. 



IV. DISCUSSION AND CONCLUSIONS 



The results discussed in this paper show that for 
SPC/E water model an analogy between the supercooled 
bulk and the confined as a function of hydration level 
of the pore is possible only to the extent that in con- 
fined water upon decreasing the hydration level a glassy 
behaviour appears already at ambient temperature while 
for bulk water supercooling is required. None-the-less the 
manner the confined liquid approaches the kinetic glass 
transition temperature appears completely different. In 
the pore at higher hydration levels, due to the hydrophilic 
surface, two quite distinct subsets of water molecules are 
detectable. Experiments in favour of two phases for con- 
fined hydrogen bonded liquids are present in literature 
p7[ . The subset that is in contact with the surface is 
at higher density with respect to the bulk and is already 
a glass with low mobility even at ambient temperature. 
The inner subset displays, like a supercooled liquid, a two 
step relaxation behaviour a la MCT. It behaves none-the- 
less differently from the supercooled bulk. In fact it turns 
out that confined and bulk relaxation times are compa- 
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